{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# 4.6 Lab: Logistic Regression, LDA, QDA, and KNN"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-07-05T12:00:37.180873Z",
     "start_time": "2023-07-05T12:00:36.707394Z"
    }
   },
   "outputs": [],
   "source": [
    "import pandas as pd\n",
    "import numpy as np\n",
    "\n",
    "import statsmodels.api as sm\n",
    "import statsmodels.formula.api as smf\n",
    "\n",
    "from sklearn.linear_model import LogisticRegression\n",
    "from sklearn.discriminant_analysis import LinearDiscriminantAnalysis as LDA\n",
    "from sklearn.discriminant_analysis import QuadraticDiscriminantAnalysis as QDA\n",
    "from sklearn.neighbors import KNeighborsClassifier\n",
    "from sklearn.preprocessing import StandardScaler\n",
    "from sklearn.metrics import accuracy_score,confusion_matrix\n",
    "\n",
    "import matplotlib.pyplot as plt\n",
    "import seaborn as sns"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4.6.1 The Stock Market Data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-07-05T12:00:37.201296Z",
     "start_time": "2023-07-05T12:00:37.182237Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(1250, 9)\n"
     ]
    },
    {
     "data": {
      "text/plain": "   Year   Lag1   Lag2   Lag3   Lag4   Lag5  Volume  Today Direction\n1  2001  0.381 -0.192 -2.624 -1.055  5.010  1.1913  0.959        Up\n2  2001  0.959  0.381 -0.192 -2.624 -1.055  1.2965  1.032        Up\n3  2001  1.032  0.959  0.381 -0.192 -2.624  1.4112 -0.623      Down\n4  2001 -0.623  1.032  0.959  0.381 -0.192  1.2760  0.614        Up\n5  2001  0.614 -0.623  1.032  0.959  0.381  1.2057  0.213        Up",
      "text/html": "<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead th {\n        text-align: right;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr style=\"text-align: right;\">\n      <th></th>\n      <th>Year</th>\n      <th>Lag1</th>\n      <th>Lag2</th>\n      <th>Lag3</th>\n      <th>Lag4</th>\n      <th>Lag5</th>\n      <th>Volume</th>\n      <th>Today</th>\n      <th>Direction</th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th>1</th>\n      <td>2001</td>\n      <td>0.381</td>\n      <td>-0.192</td>\n      <td>-2.624</td>\n      <td>-1.055</td>\n      <td>5.010</td>\n      <td>1.1913</td>\n      <td>0.959</td>\n      <td>Up</td>\n    </tr>\n    <tr>\n      <th>2</th>\n      <td>2001</td>\n      <td>0.959</td>\n      <td>0.381</td>\n      <td>-0.192</td>\n      <td>-2.624</td>\n      <td>-1.055</td>\n      <td>1.2965</td>\n      <td>1.032</td>\n      <td>Up</td>\n    </tr>\n    <tr>\n      <th>3</th>\n      <td>2001</td>\n      <td>1.032</td>\n      <td>0.959</td>\n      <td>0.381</td>\n      <td>-0.192</td>\n      <td>-2.624</td>\n      <td>1.4112</td>\n      <td>-0.623</td>\n      <td>Down</td>\n    </tr>\n    <tr>\n      <th>4</th>\n      <td>2001</td>\n      <td>-0.623</td>\n      <td>1.032</td>\n      <td>0.959</td>\n      <td>0.381</td>\n      <td>-0.192</td>\n      <td>1.2760</td>\n      <td>0.614</td>\n      <td>Up</td>\n    </tr>\n    <tr>\n      <th>5</th>\n      <td>2001</td>\n      <td>0.614</td>\n      <td>-0.623</td>\n      <td>1.032</td>\n      <td>0.959</td>\n      <td>0.381</td>\n      <td>1.2057</td>\n      <td>0.213</td>\n      <td>Up</td>\n    </tr>\n  </tbody>\n</table>\n</div>"
     },
     "execution_count": 2,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data = pd.read_csv('../data/Smarket.csv',index_col=0)\n",
    "print(data.shape)\n",
    "data.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-07-05T12:00:37.260262Z",
     "start_time": "2023-07-05T12:00:37.201942Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": "Index(['Year', 'Lag1', 'Lag2', 'Lag3', 'Lag4', 'Lag5', 'Volume', 'Today',\n       'Direction'],\n      dtype='object')"
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data.columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-07-05T12:00:37.291128Z",
     "start_time": "2023-07-05T12:00:37.206733Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": "              Year         Lag1         Lag2         Lag3         Lag4  \\\ncount  1250.000000  1250.000000  1250.000000  1250.000000  1250.000000   \nmean   2003.016000     0.003834     0.003919     0.001716     0.001636   \nstd       1.409018     1.136299     1.136280     1.138703     1.138774   \nmin    2001.000000    -4.922000    -4.922000    -4.922000    -4.922000   \n25%    2002.000000    -0.639500    -0.639500    -0.640000    -0.640000   \n50%    2003.000000     0.039000     0.039000     0.038500     0.038500   \n75%    2004.000000     0.596750     0.596750     0.596750     0.596750   \nmax    2005.000000     5.733000     5.733000     5.733000     5.733000   \n\n             Lag5       Volume        Today  \ncount  1250.00000  1250.000000  1250.000000  \nmean      0.00561     1.478305     0.003138  \nstd       1.14755     0.360357     1.136334  \nmin      -4.92200     0.356070    -4.922000  \n25%      -0.64000     1.257400    -0.639500  \n50%       0.03850     1.422950     0.038500  \n75%       0.59700     1.641675     0.596750  \nmax       5.73300     3.152470     5.733000  ",
      "text/html": "<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead th {\n        text-align: right;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr style=\"text-align: right;\">\n      <th></th>\n      <th>Year</th>\n      <th>Lag1</th>\n      <th>Lag2</th>\n      <th>Lag3</th>\n      <th>Lag4</th>\n      <th>Lag5</th>\n      <th>Volume</th>\n      <th>Today</th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th>count</th>\n      <td>1250.000000</td>\n      <td>1250.000000</td>\n      <td>1250.000000</td>\n      <td>1250.000000</td>\n      <td>1250.000000</td>\n      <td>1250.00000</td>\n      <td>1250.000000</td>\n      <td>1250.000000</td>\n    </tr>\n    <tr>\n      <th>mean</th>\n      <td>2003.016000</td>\n      <td>0.003834</td>\n      <td>0.003919</td>\n      <td>0.001716</td>\n      <td>0.001636</td>\n      <td>0.00561</td>\n      <td>1.478305</td>\n      <td>0.003138</td>\n    </tr>\n    <tr>\n      <th>std</th>\n      <td>1.409018</td>\n      <td>1.136299</td>\n      <td>1.136280</td>\n      <td>1.138703</td>\n      <td>1.138774</td>\n      <td>1.14755</td>\n      <td>0.360357</td>\n      <td>1.136334</td>\n    </tr>\n    <tr>\n      <th>min</th>\n      <td>2001.000000</td>\n      <td>-4.922000</td>\n      <td>-4.922000</td>\n      <td>-4.922000</td>\n      <td>-4.922000</td>\n      <td>-4.92200</td>\n      <td>0.356070</td>\n      <td>-4.922000</td>\n    </tr>\n    <tr>\n      <th>25%</th>\n      <td>2002.000000</td>\n      <td>-0.639500</td>\n      <td>-0.639500</td>\n      <td>-0.640000</td>\n      <td>-0.640000</td>\n      <td>-0.64000</td>\n      <td>1.257400</td>\n      <td>-0.639500</td>\n    </tr>\n    <tr>\n      <th>50%</th>\n      <td>2003.000000</td>\n      <td>0.039000</td>\n      <td>0.039000</td>\n      <td>0.038500</td>\n      <td>0.038500</td>\n      <td>0.03850</td>\n      <td>1.422950</td>\n      <td>0.038500</td>\n    </tr>\n    <tr>\n      <th>75%</th>\n      <td>2004.000000</td>\n      <td>0.596750</td>\n      <td>0.596750</td>\n      <td>0.596750</td>\n      <td>0.596750</td>\n      <td>0.59700</td>\n      <td>1.641675</td>\n      <td>0.596750</td>\n    </tr>\n    <tr>\n      <th>max</th>\n      <td>2005.000000</td>\n      <td>5.733000</td>\n      <td>5.733000</td>\n      <td>5.733000</td>\n      <td>5.733000</td>\n      <td>5.73300</td>\n      <td>3.152470</td>\n      <td>5.733000</td>\n    </tr>\n  </tbody>\n</table>\n</div>"
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "data.describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-07-05T12:00:37.291462Z",
     "start_time": "2023-07-05T12:00:37.224571Z"
    }
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/var/folders/fp/3_bvcq4d46j3zghpsc9ym8_00000gn/T/ipykernel_5619/968127052.py:2: FutureWarning: The default value of numeric_only in DataFrame.corr is deprecated. In a future version, it will default to False. Select only valid columns or specify the value of numeric_only to silence this warning.\n",
      "  corr = data.corr()\n"
     ]
    },
    {
     "data": {
      "text/plain": "            Year      Lag1      Lag2      Lag3      Lag4      Lag5    Volume  \\\nYear    1.000000  0.029700  0.030596  0.033195  0.035689  0.029788  0.539006   \nLag1    0.029700  1.000000 -0.026294 -0.010803 -0.002986 -0.005675  0.040910   \nLag2    0.030596 -0.026294  1.000000 -0.025897 -0.010854 -0.003558 -0.043383   \nLag3    0.033195 -0.010803 -0.025897  1.000000 -0.024051 -0.018808 -0.041824   \nLag4    0.035689 -0.002986 -0.010854 -0.024051  1.000000 -0.027084 -0.048414   \nLag5    0.029788 -0.005675 -0.003558 -0.018808 -0.027084  1.000000 -0.022002   \nVolume  0.539006  0.040910 -0.043383 -0.041824 -0.048414 -0.022002  1.000000   \nToday   0.030095 -0.026155 -0.010250 -0.002448 -0.006900 -0.034860  0.014592   \n\n           Today  \nYear    0.030095  \nLag1   -0.026155  \nLag2   -0.010250  \nLag3   -0.002448  \nLag4   -0.006900  \nLag5   -0.034860  \nVolume  0.014592  \nToday   1.000000  ",
      "text/html": "<div>\n<style scoped>\n    .dataframe tbody tr th:only-of-type {\n        vertical-align: middle;\n    }\n\n    .dataframe tbody tr th {\n        vertical-align: top;\n    }\n\n    .dataframe thead th {\n        text-align: right;\n    }\n</style>\n<table border=\"1\" class=\"dataframe\">\n  <thead>\n    <tr style=\"text-align: right;\">\n      <th></th>\n      <th>Year</th>\n      <th>Lag1</th>\n      <th>Lag2</th>\n      <th>Lag3</th>\n      <th>Lag4</th>\n      <th>Lag5</th>\n      <th>Volume</th>\n      <th>Today</th>\n    </tr>\n  </thead>\n  <tbody>\n    <tr>\n      <th>Year</th>\n      <td>1.000000</td>\n      <td>0.029700</td>\n      <td>0.030596</td>\n      <td>0.033195</td>\n      <td>0.035689</td>\n      <td>0.029788</td>\n      <td>0.539006</td>\n      <td>0.030095</td>\n    </tr>\n    <tr>\n      <th>Lag1</th>\n      <td>0.029700</td>\n      <td>1.000000</td>\n      <td>-0.026294</td>\n      <td>-0.010803</td>\n      <td>-0.002986</td>\n      <td>-0.005675</td>\n      <td>0.040910</td>\n      <td>-0.026155</td>\n    </tr>\n    <tr>\n      <th>Lag2</th>\n      <td>0.030596</td>\n      <td>-0.026294</td>\n      <td>1.000000</td>\n      <td>-0.025897</td>\n      <td>-0.010854</td>\n      <td>-0.003558</td>\n      <td>-0.043383</td>\n      <td>-0.010250</td>\n    </tr>\n    <tr>\n      <th>Lag3</th>\n      <td>0.033195</td>\n      <td>-0.010803</td>\n      <td>-0.025897</td>\n      <td>1.000000</td>\n      <td>-0.024051</td>\n      <td>-0.018808</td>\n      <td>-0.041824</td>\n      <td>-0.002448</td>\n    </tr>\n    <tr>\n      <th>Lag4</th>\n      <td>0.035689</td>\n      <td>-0.002986</td>\n      <td>-0.010854</td>\n      <td>-0.024051</td>\n      <td>1.000000</td>\n      <td>-0.027084</td>\n      <td>-0.048414</td>\n      <td>-0.006900</td>\n    </tr>\n    <tr>\n      <th>Lag5</th>\n      <td>0.029788</td>\n      <td>-0.005675</td>\n      <td>-0.003558</td>\n      <td>-0.018808</td>\n      <td>-0.027084</td>\n      <td>1.000000</td>\n      <td>-0.022002</td>\n      <td>-0.034860</td>\n    </tr>\n    <tr>\n      <th>Volume</th>\n      <td>0.539006</td>\n      <td>0.040910</td>\n      <td>-0.043383</td>\n      <td>-0.041824</td>\n      <td>-0.048414</td>\n      <td>-0.022002</td>\n      <td>1.000000</td>\n      <td>0.014592</td>\n    </tr>\n    <tr>\n      <th>Today</th>\n      <td>0.030095</td>\n      <td>-0.026155</td>\n      <td>-0.010250</td>\n      <td>-0.002448</td>\n      <td>-0.006900</td>\n      <td>-0.034860</td>\n      <td>0.014592</td>\n      <td>1.000000</td>\n    </tr>\n  </tbody>\n</table>\n</div>"
     },
     "execution_count": 5,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "#produces a correlation matrix for all the numrice columns \n",
    "corr = data.corr()\n",
    "corr"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-07-05T12:00:37.605629Z",
     "start_time": "2023-07-05T12:00:37.235566Z"
    }
   },
   "outputs": [
    {
     "data": {
      "text/plain": "<AxesSubplot: >"
     },
     "execution_count": 6,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": "<Figure size 800x800 with 2 Axes>",
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAnQAAAKWCAYAAADA5QRGAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjYuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8o6BhiAAAACXBIWXMAAA9hAAAPYQGoP6dpAABJl0lEQVR4nO3deXwV9b3/8fdJyA7IKgEMRTaXogYTQBSkBARJJSxGUBCFuiAQFqFWjOwSlCsqIsVGW4sKBSwqRgGtCnrhRkgEtZrKEm6BaEokCSFkgWzz+8MfuR4mkDMhh8mQ17OP8yhnzpwzn0kwfPL+zvc7LsMwDAEAAMCxfOwuAAAAABeGhg4AAMDhaOgAAAAcjoYOAADA4WjoAAAAHI6GDgAAwOFo6AAAAByOhg4AAMDhaOgAAAAcjoYOAADAi3Jzc3Xbbbdp165d59zn888/15AhQxQeHq7Bgwdr27Ztlo5BQwcAAOAlu3fv1qhRo3TkyJFz7nPo0CFNmTJF06ZN05dffqkpU6Zo+vTpysrK8vg4NHQAAABe8O677+r3v/+9Hn300Wr3i4yM1IABA9SgQQNFR0ere/fuWr9+vcfHoqEDAADwUElJiQoKCtweJSUlVe7bu3dvffzxx4qOjj7vZ6anp6tLly5u2zp16qS9e/d6XFcDj/f0ktLs/7W7BFt07DLU7hIAr3O5XHaXYAvDMOwuwRafh4XaXYIt+mYctbsEWxzJ/da2Y9vZOySu3aQVK1a4bYuLi9OUKVNM+7Zs2dKjzywsLFRQUJDbtsDAQBUVFXlcl+0NHQAAgFNMmDBB48ePd9vm7+9/QZ8ZFBSkU6dOuW07deqUQkJCPP4MGjoAAAAP+fv7X3ADd7YuXbooLS3NbVt6erq6du3q8WdwDR0AAHCWinL7Hl4QExOjlJQUbd68WWVlZdq8ebNSUlI0dKjnl2fR0AEAAFxk3bp1U1JSkiSpY8eO+uMf/6jExER1795dK1eu1EsvvaQrr7zS489jyBUAADiLUWF3BZbt27fP7flXX33l9rxPnz7q06dPjT+fhA4AAMDhSOgAAICzVDgvofM2EjoAAACHo6EDAABwOIZcAQCAoxgOnBThbSR0AAAADkdCBwAAnIVJESYkdAAAAA5HQwcAAOBwDLkCAABnYVKECQkdAACAw5HQAQAAZ6kot7uCOoeEDgAAwOFo6AAAAByOIVcAAOAsTIowIaEDAABwOBI6AADgLNwpwoSEDgAAwOFI6AAAgKMYXENnQkIHAADgcJYauoyMDG/VAQAAgBqy1NCNGjVKBQUF3qoFAACgehUV9j3qKEsNXZMmTZSVleWtWgAAAFADliZFdO7cWSNHjlR4eLguv/xyt9eefvrpWi0MAACgSkyKMLHU0AUHB2vgwIHeqgUAAAA1YKmhI4UDAACoeyw1dCUlJXr//feVlZWliv9/YWBpaan279+vl19+2SsFAgAAuKkot7uCOsdSQxcfH6/t27eradOmKi0tVXBwsA4cOKBhw4Z5qTwAAABUx1JDt337dq1du1a5ublau3atnnvuOb322mv65z//6a36AAAA3DEpwsTSsiUVFRXq0KGDOnTooO+//16SNGbMGH355ZdeKQ4AAADVs5TQhYaGKiMjQ2FhYcrJyVFRUZF8fHxUWFjorfoAAADc1eEFfu1iqaEbMmSIRo8erQ0bNug3v/mNJk6cqICAAHXt2tVb9QEAAKAalhq6hx9+WGFhYWrUqJHmzJmjZ599VgUFBZozZ4636gMAAEA1LDV0kjR48GBJUm5urhYsWFDrBQEAAJwXkyJMLE2KKC0t1QsvvKCIiAhFRUUpIyNDd955p3766Sdv1QcAAIBqWGroVqxYoZ07d+rFF1+Un5+fmjdvrtDQUCUkJHirPgAAAHcVFfY96ihLQ67vv/++1q5dq1atWsnlcik4OFhPP/20brvtNm/VBwAAgGpYSuiKiorUrFkzSZJhGJKkwMBA+fhY+hgAAADUIkudWHh4uFasWCFJcrlckqQ333xT1113Xe1XBgAAUAXDKLftUVd51NA98sgjOnnypJ588km9//77uvXWW1VYWKjo6Gi98cYbmjVrlrfrBAAAwDl4dA1dbm6uhg0bpuXLl2vTpk3atm2bMjMzFRoaqt/85jdq2LCht+sEAAD4GcuWmHjU0K1du1bLly/Xvffeq8cff1x33323t+sCAACAhzxq6Hx9ffXoo4/qlltu0RNPPKHdu3dr2rRpbpMh2rRp47UiAQAAKtXh5UPsYmnZkh49euiPf/yj7r77bn3wwQeSfp7t6nK59P3333ulQAAAAJyfpYZuzZo1Wrp0qQYNGqS4uDiWKwEAAKgDPJ4UER8fr5SUFM2bN0/Dhg3zclkAAADnwKQIE48aupiYGLVs2VJvv/22rrzySm/XBAAAAAs8auhuu+02PfHEE/L39/d2PQAAAOdXUXcX+LWLRw3dvHnzvF0HAAAAaohZDQAAAA5naZYrAACA7ZgUYUJCBwAA4HAkdAAAwFm4U4QJCR0AAIDDkdABAABn4Ro6ExI6AAAAh6OhAwAAcDiGXAEAgLMwKcKEhA4AAMDhSOgAAICzkNCZkNABAAA4HA0dAACAwzHkCgAAHMUwyu0uoc4hoQMAAHA4EjoAAOAsTIowIaEDAABwOBI6AADgLNzL1YSEDgAAwOFo6AAAAByOIVcAAOAsTIowIaEDAABwONsTuo5dhtpdgi0O7n/P7hJs0Tisn90l2MLf1/b/1Gzh66qfvzMahmF3CbaI+uEnu0uwxbHifLtLqH+YFGFSP3/aAgAAXEJo6AAAAByufo4DAQAA52JShAkJHQAAgMOR0AEAAGdhUoQJCR0AAIDDkdABAABn4Ro6ExI6AAAAh6OhAwAAcDiGXAEAgLMw5GpCQgcAAOBwJHQAAMBZWLbEhIQOAADA4WjoAAAAHI4hVwAA4CxMijAhoQMAAHA4EjoAAOAsTIowIaEDAABwOBo6AAAAh2PIFQAAOAuTIkxI6AAAAByOhA4AADgLkyJMSOgAAAAcjoQOAAA4C9fQmZDQAQAAOBwNHQAAgMN5POS6cePGavcZNmzYBZQCAADgAYZcTTxu6NatW6dvvvlGrVu3rvJ1l8tFQwcAAGADjxu6v/71r7r33ns1YsQIjRkzxps1AQAAnJth2F1BnePxNXRBQUF65pln9PLLL6ukpMSbNQEAAMACS5MiOnfurKVLl6q4uNhb9QAAAMAiy+vQ3XTTTd6oAwAAwDNMijCx3NClpqZWud3Pz0/NmjVTu3btLrgoAAAAeM5yQzdr1ixlZmbKx8dHTZs21fHjx1VRUSEfHx+Vl5erQ4cOSkxMVFhYmDfqBQAA9R0JnYnlhYVjYmIUExOjlJQU7dixQ6mpqYqNjVVcXJx2796t3r17KyEhwRu1AgAAoAqWG7qNGzdq/vz5CgkJkSQFBwcrPj5e69evV0hIiGbOnKk9e/bUeqEAAACSJKPCvkcdZbmhKyoqUn5+vtu2kydPqqCgoPK5y+W68MoAAAAcLCcnR5MmTVJkZKR69uyphIQElZWVVbnv66+/rqioKN14440aMmSIPvroI0vHstzQ3X777Zo8ebKSk5N16NAhJScna+rUqRo4cKAKCgo0b948RUZGWv1YAACAS8r06dMVHBys7du3a8OGDfriiy+0atUq036ff/65EhMT9ec//1l79uxRXFycpk+frh9++MHjY1meFBEfH6+EhARNnjxZxcXFCgwMVGxsrGbOnKm0tDTl5+dr/vz5Vj8WAADAMw6YFHH48GGlpKTov//7vxUUFKSwsDBNmjRJzz77rB588EG3ff/3f/9XhmFUPnx9feXn56cGDTxv0yw3dAEBAVq4cKHmzp2rvLw8NW/evHKINTIyknQOAABcskpKSkx3zPL395e/v7/btgMHDqhJkyZq1apV5baOHTsqMzNT+fn5aty4ceX23/72t3rnnXcUHR0tX19fuVwuPfvsswoNDfW4LssNnSTt3LlTWVlZMv7/vdRKS0u1b98+zZ49uyYfBwAA4Dkb7+WamJioFStWuG2Li4vTlClT3LYVFhYqKCjIbduZ50VFRW4NXWlpqa6++molJCTo6quv1vvvv68nn3xSHTt21FVXXeVRXZYbukWLFmndunWVs1zLy8tVWFioPn36WP0oAAAAR5kwYYLGjx/vtu3sdE76eRWQs2+Veub5mR7qjKeeeko33nijrr/+eknSnXfeqQ8++EDvvvuuZs2a5VFdlhu6LVu2aPXq1SouLlZSUpIWL16sJUuWqKioyOpHAQAAOEpVw6tV6dy5s/Ly8pSdna0WLVpIkg4ePKjQ0FA1atTIbd/MzEx17drVbVuDBg3k5+fncV2WZ7kWFxcrPDxcnTp1Ulpamlwul+Li4vTZZ59Z/SgAAADrKirse3ioffv2ioiI0OLFi1VQUKCMjAytXLlSsbGxpn2joqK0evVqpaWlqaKiQh9++KF27dql6Ohoj49nOaELDQ1VTk6OWrZsqaNHj6q0tFSBgYFu69ABAADUd8uXL9fChQvVv39/+fj4aNiwYZo0aZIkqVu3blqwYIFiYmIUFxcnX19fTZkyRSdOnNCvfvUr/fGPf9Q111zj8bFchmHtysIlS5Zox44dev311zVnzhwFBwcrICBA//rXv/TOO+9YO1NJ7ZpdZ/k9l4KD+9+zuwRbNA7rZ3cJtvD3rdH8I8fzdVkeBLgkWPyxesloEhBS/U6XoKNFeXaXYIvi4sP2Hfsvv7ft2EEPLLXt2Odj+aftjBkzNHToUPn5+VUuXZKenq6nnnrKG/UBAACgGpZjAz8/v8oF8Ro1aqRXX31V5eXlOnLkSK0XBwAAYFKH76lql1oZD8nOzrZ04R4AAABqT61d4FJfrxkBAACwW61dqX3m9l8AAADeZFQQIp2tfk5BAwAAuIR4nNClpqae87Xc3NxaKQYAAKBaFhb4rS88bujGjh173tcZcgUAALCHxw3d3r17vVkHAAAAaqh+Ll8PAACci3XoTJgUAQAA4HAkdAAAwFlYtsSEhA4AAMDhSOgAAICzsGyJCQkdAACAw9HQAQAAOBxDrgAAwFkYcjUhoQMAAHA4EjoAAOAsBsuWnI2EDgAAwOFo6AAAAByOIVcAAOAsTIowIaEDAABwOBI6AADgLNzL1YSEDgAAwOFI6AAAgLMYXEN3NhI6AAAAh6OhAwAAcDiGXAEAgLMwKcKEhA4AAMDhSOhs0jisn90l2CI/Y5vdJdiivn6//X3r548YH7nsLsEWP5zMtrsEWwT7B9pdQr1jsLCwCQkdAACAw9HQAQAAOFz9HA8BAADOxaQIExI6AAAAhyOhAwAAzsKdIkxI6AAAAByOhA4AADgL19CZkNABAAA4HA0dAACAwzHkCgAAnIU7RZiQ0AEAADgcCR0AAHAWJkWYkNABAAA4HA0dAACAwzHkCgAAnIU7RZiQ0AEAADgcCR0AAHAWJkWYkNABAAA4HA0dAACAwzHkCgAAHMXgThEmJHQAAAAOR0IHAACchUkRJiR0AAAADkdCBwAAnIWEzsRyQnfy5Mkqtx89evSCiwEAAIB1Hjd0P/zwg4YMGaIePXpowIAB2rlzp9vr0dHRtV4cAAAAqudxQ/fMM8/o+uuv18aNGzVkyBBNnDhRqampla8bBvEnAAC4CIwK+x51lMfX0O3evVtbt25VUFCQrrrqKl1++eWaMmWK3nnnHbVp00Yul8ubdQIAAOAcLF1D5+Pzf7vfc889Gjx4sKZMmaKSkhISOgAAcHFUGPY96iiPG7rIyEglJCQoOzu7clt8fLz8/Pw0efJkGjoAAACbeNzQzZo1S//85z/1xBNPVG7z8/PTyy+/rJycHJ0+fdorBQIAAOD8PL6Grm3bttq4caPy8/Pdtjdt2lRvvfWWtm3bVuvFAQAAnM2ow0OfdrG8sPC+ffuq3N6yZUsdOXJE7dq1u+CiAAAA4DnLDd2sWbOUmZkpHx8fNW3aVMePH1dFRYV8fHxUXl6uDh06KDExUWFhYd6oFwAA1HckdCaW7xQRExOjmJgYpaSkaMeOHUpNTVVsbKzi4uK0e/du9e7dWwkJCd6oFQAAAFWwnNBt3LhRmzdvVlBQkCQpODhY8fHxGjx4sCZOnKiZM2eqT58+tV4oAACAJKmi7i7waxfLCV1RUZFpYsTJkydVUFBQ+ZxFhgEAAC4eyw3d7bffrsmTJys5OVmHDh1ScnKypk6dqoEDB6qgoEDz5s1TZGSkN2oFAABAFSwPucbHxyshIUGTJ09WcXGxAgMDFRsbq5kzZyotLU35+fmaP3++F0oFAAAQkyKq4DJqeIuHsrIy5eXlqXnz5hc0xNqu2XU1fq+THSvOr36nS1B+Rv1cr7BxWD+7S7CFv6/l3xkvCT6qn5edFJQU212CLYL9A+0uwRYnCg7aduyTkwbbduxGK7fYduzzqdFP2507dyorK6vydl+lpaXat2+fZs+eXavFAQAAmJDQmVhu6BYtWqR169YpJCREklReXq7CwkJmtgIAANjEckO3ZcsWrV69WsXFxUpKStLixYu1ZMkSFRUVeaM+AAAAVMNyQ1dcXKzw8HAdO3ZMaWlpcrlciouLU3R0tDfqAwAAcFPDy/8vaZaXLQkNDVVOTo5atmypo0ePqrS0VIGBgW7r0AEAAODisZzQ9e3bV+PGjdPrr7+u7t27Kz4+XgEBAWrfvr0XygMAADgLkyJMLCd0M2bM0NChQ+Xn56e5c+cqLy9P6enpeuqpp7xRHwAAAKphOaHz8/PTgw8+KElq1KiRXn31VZWXl+vIkSO1XhwAAIAJCZ2J5YSuKtnZ2UyKAAAAsEmtNHQSM04AAADsUmv35bmQ238BAAB4ymDI1aTWEjoAAADYw+OELjU19Zyv5ebm1koxAAAA1SKhM/G4oRs7dux5X2fIFQAAwB4eN3R79+71Zh0AAACooVqbFAEAAHBRVNhdQN3DpAgAAACHI6EDAACOwrIlZiR0AAAADkdCBwAAnIWEzoSEDgAAwOFo6AAAAByOIVcAAOAsLFtiQkIHAADgcCR0AADAUVi2xIyEDgAAwOFo6AAAAByOIVcAAOAsTIowIaEDAABwOBI6AADgKEyKMCOhAwAAcDgSOgAA4CxcQ2dCQgcAAOBwNHQAAAAOx5ArAABwFIMhVxMSOgAAAIcjobOJv2/9/NI3Dutndwm2yM/YZncJtqiv329fF78r1yd8v23gkIQuJydHc+bMUUpKinx9fRUTE6PHH39cDRqYe4CUlBQ9++yzSk9PV+PGjTV69GhNmDDB42PxtxAAAMALpk+fruDgYG3fvl0bNmzQF198oVWrVpn2O3jwoB5++GGNHj1ae/bsUWJiol577TV9+OGHHh+Lhg4AAKCWHT58WCkpKXrssccUFBSksLAwTZo0SWvWrDHt+7e//U39+/fX8OHD5XK5dPXVV2vdunWKiIjw+Hg0dAAAwFGMCvseJSUlKigocHuUlJSYajxw4ICaNGmiVq1aVW7r2LGjMjMzlZ+f77bvP//5T11xxRWaMWOGevbsqcGDByslJUUtW7b0+GtCQwcAAOChxMRERUREuD0SExNN+xUWFiooKMht25nnRUVFbttPnDihN954QzExMfqf//kfLVy4UEuWLLE05Fo/r8wHAADOZeOkiAkTJmj8+PFu2/z9/U37BQcHq7i42G3bmechISGm9/fv31+/+c1vJEndu3fX0KFDtWXLFt1+++0e1UVDBwAA4CF/f/8qG7izde7cWXl5ecrOzlaLFi0k/Tz5ITQ0VI0aNXLbt2PHjqZh2/LychmG4XFdDLkCAABHsfMaOk+1b99eERERWrx4sQoKCpSRkaGVK1cqNjbWtO/dd9+tTz/9VO+9954Mw1Bqaqref/99DR061OPj0dABAAB4wfLly1VWVqb+/ftr5MiR6tOnjyZNmiRJ6tatm5KSkiRJvXr10sqVK/XGG28oIiJCTzzxhB5//HH179/f42O5DCt5nhe0a3adnYe3zYmSoup3ugSVlJfZXYItWFi4fqmvC82eLjPP9KsPGgUE212CLXJPHrDt2Mdu62vbsVt+/Lltxz4frqEDAACOwr1czernr5EAAACXEBI6AADgKCR0ZiR0AAAADkdDBwAA4HAMuQIAAGcxXHZXUOeQ0AEAADgcCR0AAHAUJkWYkdABAAA4HA0dAACAwzHkCgAAHMWoYFLE2UjoAAAAHI6EDgAAOAqTIsxI6AAAAByOhA4AADiKwcLCJjVO6AzD0N69e/Wf//ynNusBAACARR4ndHl5eVqyZImaNWumBx54QPfff78OHDggl8ulvn37aunSpWrYsKE3awUAAEAVPE7onn76aR09elS7d+/W/fffr86dOys5OVlbt26VJC1btsxbNQIAAFQyKux71FUeJ3T//d//rY8//ljHjx/XwIED9be//U2NGjWSJC1evFjDhw/X7NmzvVYoAAAAquZxQ1daWqqQkBA1bNhQ4eHhCggIqHwtMDBQp06d8kqBAAAAv8TCwmYeD7n++te/VmJioiRp7dq18vf3lyRlZ2drzpw56tmzp3cqBAAAwHl53NDNmjVLa9asMSVxo0aN0qFDh/Tkk0/WenEAAAConsdDrtdcc422bt0qPz8/t+1vvPGG2rZtW+uFAQAAVMUw7K6g7rG0sLCfn59SU1NN2zMzM+Xn56dmzZqpXbt2tVYcAAAAqmf5ThGzZs1SZmamfHx81LRpUx0/flwVFRXy8fFReXm5OnTooMTERIWFhXmjXgAAUM8xKcLM8p0iYmJiFBMTo5SUFO3YsUOpqamKjY1VXFycdu/erd69eyshIcEbtQIAAKAKlhO6jRs3avPmzQoKCpIkBQcHKz4+XoMHD9bEiRM1c+ZM9enTp9YLBQAAkEjoqmI5oSsqKlJ+fr7btpMnT6qgoKDyucvFFxoAAOBisdzQ3X777Zo8ebKSk5N16NAhJScna+rUqRo4cKAKCgo0b948RUZGeqNWAAAAVMHykGt8fLwSEhI0efJkFRcXKzAwULGxsZo5c6bS0tKUn5+v+fPne6FUAAAAli2pisswavZlKSsrU15enpo3b35BQ6ztml1X4/c62YmSIrtLsEVJeZndJdgiP2Ob3SXYonFYP7tLsIWvy/LgxyXhdFmJ3SXYolFAsN0l2CL35AHbjv3vG26z7dhXfvOxbcc+H8sJnSTt3LlTWVlZOtMLlpaWat++fZo9e3atFgcAAHA2JkWYWW7oFi1apHXr1ikkJESSVF5ersLCQma2AgAA2MRyQ7dlyxatXr1axcXFSkpK0uLFi7VkyRIVFdXPIUQAAAC7WW7oiouLFR4ermPHjiktLU0ul0txcXGKjo72Rn0AAABuDIMh17NZvnI3NDRUOTk5atmypY4eParS0lIFBga6rUMHAACAi8dyQte3b1+NGzdOr7/+urp37674+HgFBASoffv2XigPAADAnVFhdwV1j+WEbsaMGRo6dKj8/Pw0d+5c5eXlKT09XU899ZQ36gMAAEA1LCd0fn5+evDBByVJjRo10quvvqry8nIdOXKk1osDAAA4WwXX0JnUyuqX2dnZTIoAAACwSa0tZ17DG04AAADgAtXoThFVuZDbfwEAAHiKZUvM6ucNBwEAAC4hHid0qamp53wtNze3VooBAACoDvdyNfO4oRs7dux5X2fIFQAAwB4eN3R79+71Zh0AAACooVqbFAEAAHAxsLCGGZMiAAAAHI6EDgAAOAqTIsxI6AAAAByOhA4AADgK93I1I6EDAABwOBo6AAAAh2PIFQAAOAr3cjUjoQMAAHA4EjoAAOAoLCxsRkIHAADgcDR0AAAADseQKwAAcBTWoTMjoQMAAHA4EjoAAOAoLFtiRkIHAADgcCR0AADAUVi2xIyEDgAAwOFo6AAAAByOIVcAAOAoLFtiRkIHAADgcLYndC5X/eyyfV31s5f297X9r5wtGof1s7sEW+RnbLO7BFvU1+93fb1O3eAK/YuOZUvM6mdXAQAAcAmhoQMAAHC4+jn+BQAAHItJEWYkdAAAAA5HQgcAAByFaShmJHQAAAAOR0IHAAAchWvozEjoAAAAHI6GDgAAwOEYcgUAAI7CnSLMSOgAAAAcjoQOAAA4SoXdBdRBJHQAAAAOR0MHAADgcAy5AgAARzHEpIizkdABAAA4HAkdAABwlApu5mpCQgcAAOBwNHQAAAAOx5ArAABwlAomRZiQ0AEAADgcCR0AAHAUli0xI6EDAABwOBI6AADgKNzL1eyCE7qjR4/WRh0AAACooQtu6GJiYmqjDgAAANSQx0Ou9913X5XbCwsLK1974403aqcqAACAc2BShJnHDV2rVq20adMmjRw5Ui1atKjc/s0336hHjx5eKQ4AAADV87ihe/bZZ9WzZ0+9/PLLmj9/vvr06SPp51QuLi7OawUCAAD8EpMizCzNco2NjVV4eLgeffRRffHFF5oxY4a36gIAAICHLE+K6NSpk/7+97/rxIkTGjVqlMrLy71RFwAAADxUo3XoAgMDlZCQoKSkJL333nu1XRMAAMA5MeRqZrmhS01Nrfxz69at9cgjjyg1NVV+fn5q1qyZ2rVrV6sFAgAA4PwsN3SzZs1SZmamfHx81LRpUx0/flwVFRXy8fFReXm5OnTooMTERIWFhXmjXgAAUM+xbImZ5WvoYmJiFBMTo5SUFO3YsUOpqamKjY1VXFycdu/erd69eyshIcEbtQIAAKAKlhu6jRs3av78+QoJCZEkBQcHKz4+XuvXr1dISIhmzpypPXv21HqhAAAAklThsu9RV1lu6IqKipSfn++27eTJkyooKKh87nLV4TMGAAC4xFhu6G6//XZNnjxZycnJOnTokJKTkzV16lQNHDhQBQUFmjdvniIjI71RKwAAAKpgeVJEfHy8EhISNHnyZBUXFyswMFCxsbGaOXOm0tLSlJ+fr/nz53uhVAAAAKmCSREmLsMwjJq8saysTHl5eWrevPkFDbH+qvn1NX6vk50sKba7BFuUG/Vz9aCS8jK7S7BFfsY2u0uwReOwfnaXYIvSevr3vJF/kN0l2OJ4Qbptx34vdLRtxx569G+2Hft8arSw8M6dO5WVlaUzvWBpaan27dun2bNn12pxAAAAZ6tREmWDnJwczZkzRykpKfL19VVMTIwef/xxNWhw7vZr//79uuuuu/TKK6+oZ8+eHh/LckO3aNEirVu3rnKWa3l5uQoLC9WnTx+rHwUAAHDJmj59ulq1aqXt27crOztbEydO1KpVq/Tggw9WuX9xcbFmzpypU6dOWT6W5UkRW7Zs0erVq7Vs2TJFRUUpNTVV999/v0JDQy0fHAAA4FJ0+PBhpaSk6LHHHlNQUJDCwsI0adIkrVmz5pzvWbBggQYMGFCj41lu6IqLixUeHq5OnTopLS1NLpdLcXFx+uyzz2pUAAAAgBUVNj5KSkpUUFDg9igpKTHVeODAATVp0kStWrWq3NaxY0dlZmaaln+Tfl7n9/Dhw4qLi6vR18RyQxcaGqqcnBy1bNlSR48eVWlpqQIDA93WoQMAALgUJSYmKiIiwu2RmJho2q+wsFBBQe4TZs48Lyoqctt+8OBBvfDCC3ruuefk6+tbo7osX0PXt29fjRs3Tq+//rq6d++u+Ph4BQQEqH379jUqAAAAwIoKG29gMGHCBI0fP95tm7+/v2m/4OBgFRe7r2hx5vmZeQiSdPr0aT366KOKj49XmzZtalyX5YRuxowZGjp0qPz8/DR37lzl5eUpPT1dTz31VI2LAAAAcAJ/f381bNjQ7VFVQ9e5c2fl5eUpOzu7ctvBgwcVGhqqRo0aVW779ttvdejQIT355JOKjIysvDnDI488Ymld3xqvQ/dL5eXlOnLkiK688krL72UduvqFdejqF9ahq19Yh65+sXMdur+3HmPbse/6z7knNZxt9OjRCg0N1cKFC3X8+HFNnDhRgwYN0pQpU877vquuukpvvPGGpWVLLCd0VcnOzlZ0dHRtfBQAAMAlYfny5SorK1P//v01cuRI9enTR5MmTZIkdevWTUlJSbV2rBotLFyVWgj6AAAALhktWrTQ8uXLq3ztq6++Ouf79u3bZ/lYtdbQXcjtvwAAADxVPy/eOb9aGXIFAACAfTxO6FJTU8/5Wm5ubq0UAwAAUJ0KBgVNPG7oxo4de97XGXIFAACwh8cN3d69e71ZBwAAAGqo1iZFAAAAXAwVYlTwbEyKAAAAcDgSOgAA4CisfGtGQgcAAOBwJHQAAMBRWLbEjIQOAADA4WjoAAAAHI4hVwAA4Cjcy9WMhA4AAMDhSOgAAICjsGyJGQkdAACAw9HQAQAAOBxDrgAAwFFYh86MhA4AAMDhSOgAAICjsGyJGQkdAACAw5HQAQAARyGhMyOhAwAAcDgaOgAAAIdjyBUAADiKwbIlJiR0AAAADmd7QmcY9fOObPX1vH1UP3+t8nXVz9+dGof1s7sEW+RnbLO7BFs0aRdldwm2KCo7bXcJ9Q6TIszq578yAAAAlxAaOgAAAIezfcgVAADACoZczUjoAAAAHI6EDgAAOEr9nFZ4fiR0AAAADkdCBwAAHKWifq6AdV4kdAAAAA5HQwcAAOBwDLkCAABHYdkSMxI6AAAAhyOhAwAAjkJCZ0ZCBwAA4HA0dAAAAA7HkCsAAHAU7hRhRkIHAADgcCR0AADAUbhThBkJHQAAgMPR0AEAADgcQ64AAMBRWIfOjIQOAADA4UjoAACAo7BsiRkJHQAAgMOR0AEAAEepIKMzqXFCZxiG0tPTlZGRUZv1AAAAwCKPG7qxY8dW/jkrK0vDhg3THXfcoYEDB+qBBx5Qfn6+VwoEAADA+Xnc0KWlpVX++ZlnnlHHjh2VnJysbdu2KSAgQIsXL/ZKgQAAAL9UYeOjrvL4GjrD+L/x6tTUVH3wwQdq0qSJJGnRokUaPHhwrRcHAACA6nnc0Llc/3fjtJCQEDVo8H9vbdiwoVvDBwAA4C10HGYeD7kWFxfr/vvv17PPPqvWrVtr3bp1lduXLl2qrl27eq1IAAAAnJvHDd369es1cOBA5eTk6NixY0pOTpYkvfjii9q8ebNmzZrltSIBAABwbh4PuV5//fW6/vrrK59XVPx8aeC4ceP06KOPKiAgoParAwAAOEtdnpxgF8sLC6emppq2ZWRkyM/PT82aNVO7du1qpTAAAAB4xnJDN2vWLGVmZsrHx0dNmzbV8ePHVVFRIR8fH5WXl6tDhw5KTExUWFiYN+oFAAD1XIWr+n3qG8t3ioiJiVFMTIxSUlK0Y8cOpaamKjY2VnFxcdq9e7d69+6thIQEb9QKAACAKlhu6DZu3Kj58+crJCREkhQcHKz4+HitX79eISEhmjlzpvbs2VPrhQIAAEg/38vVrkddZbmhKyoqMt3m6+TJkyooKKh8/ss16wAAAOBdlhu622+/XZMnT1ZycrIOHTqk5ORkTZ06VQMHDlRBQYHmzZunyMhIb9QKAACAKlieFBEfH6+EhARNnjxZxcXFCgwMVGxsrGbOnKm0tDTl5+dr/vz5XigVAACAO0VUxWXU8J5dZWVlysvLU/PmzS9oiLVds+tq/F4nO1lSbHcJuIhKK8rtLsEW5Ub9XC0qP2Ob3SXYokm7KLtLsEVZPf3v+/SpDNuO/WT70bYdO+HQ32w79vlYTugkaefOncrKyqq8f2tpaan27dun2bNn12pxAAAAZ6ufvyqen+WGbtGiRVq3bl3lLNfy8nIVFhaqT58+tV4cAAAAqme5oduyZYtWr16t4uJiJSUlafHixVqyZImKioq8UR8AAACqYbmhKy4uVnh4uI4dO6a0tDS5XC7FxcUpOjraG/UBAAC4qcvrwdnF8rIloaGhysnJUcuWLXX06FGVlpYqMDDQbR06AAAAXDyWE7q+fftq3Lhxev3119W9e3fFx8crICBA7du390J5AAAA7sjnzCwndDNmzNDQoUPl5+enuXPnKi8vT+np6Xrqqae8UR8AAACqYTmh8/Pz04MPPihJatSokV599VWVl5fryJEjtV4cAADA2Vi2xMxyQleV7OxsJkUAAADYpFYaOkmq4Q0nAAAAcIFqdKeIqlzI7b8AAAA8xbIlZrWW0AEAAMAeHid0qamp53wtNze3VooBAACoDvmcmccN3dixY8/7OkOuAAAA9vC4odu7d6836wAAAEAN1dqkCAAAgIuBdejMmBQBAADgcCR0AADAUQymRZiQ0AEAADgcCR0AAHAUrqEzI6EDAABwOBo6AAAAh2PIFQAAOAr3cjUjoQMAAHA4EjoAAOAo5HNmJHQAAAAOR0MHAADgcAy5AgAAR2FShBkJHQAAgMOR0AEAAEfhThFmJHQAAAAOR0IHAAAcxeAaOhMSOgAAAIejoQMAAHA4hlwBAICjMCnCjIQOAADA4WxP6D4PC7W7BFtE/fCT3SXY4oeT2XaXgIuovl623KRdlN0l2CLvyFa7S7BFWKff2l1CvcOkCDMSOgAAAIejoQMAAHA424dcAQAArGBShBkJHQAAgBfk5ORo0qRJioyMVM+ePZWQkKCysrIq9127dq0GDRqkbt26adCgQVqzZo2lY5HQAQAAR6kwnDEpYvr06WrVqpW2b9+u7OxsTZw4UatWrdKDDz7ott8nn3yi559/Xq+++qpuuOEGff3113r44YfVokULDRo0yKNjkdABAADUssOHDyslJUWPPfaYgoKCFBYWpkmTJlWZvGVlZemhhx5SeHi4XC6XunXrpp49eyo1NdXj45HQAQAAR7EznyspKVFJSYnbNn9/f/n7+7ttO3DggJo0aaJWrVpVbuvYsaMyMzOVn5+vxo0bV24fM2aM23tzcnKUmpqqJ554wuO6SOgAAAA8lJiYqIiICLdHYmKiab/CwkIFBQW5bTvzvKio6Jyff+zYMT300EPq2rWr7rjjDo/rIqEDAADw0IQJEzR+/Hi3bWenc5IUHBys4uJit21nnoeEhFT52V9//bWmTZumyMhIPf3002rQwPM2jYYOAAA4SoWNg65VDa9WpXPnzsrLy1N2drZatGghSTp48KBCQ0PVqFEj0/4bNmzQokWLNHXqVP3ud7+zXBdDrgAAALWsffv2ioiI0OLFi1VQUKCMjAytXLlSsbGxpn0/+ugjzZ8/Xy+99FKNmjmJhg4AADiMYeP/rFi+fLnKysrUv39/jRw5Un369NGkSZMkSd26dVNSUpIkacWKFSovL9fUqVPVrVu3ysfcuXM9PhZDrgAAAF7QokULLV++vMrXvvrqq8o/v//++xd8LBI6AAAAhyOhAwAAjsK9XM1I6AAAAByOhA4AADiKncuW1FWWE7qDBw9q0aJFiouL0/Hjx7V69Wpv1AUAAAAPWWro/ud//kcjR47U8ePHlZycrFOnTumPf/yjXnnlFW/VBwAAgGpYauief/55Pf/883ruuefk6+ur1q1b65VXXtH69eu9VR8AAIAbp6xDdzFZaugOHz6sW2+9VZLkcrkkSdddd51OnDhR+5UBAADAI5YaujZt2mjPnj1u27799lu1bt26VosCAAA4lwobH3WVpVmuEyZM0MSJE3XPPfeotLRUr776qt58803NmDHDW/UBAACgGpYaut/+9rdq2LCh1qxZozZt2mjnzp168sknNWjQIG/VBwAA4MYw6u61bHaxvA5d37591bdvX2/UAgAAgBqw1NBlZGToT3/6k3788UdVVLiPJL/xxhu1WhgAAAA8Y6mhmzFjhvz8/HTTTTfJx4e7hgEAgIuPO0WYWWro0tPT9cUXXygwMNBb9QAAAMAiSzHb1VdfraNHj3qrFgAAgGqxbImZpYRu9uzZGjdunAYOHKjGjRu7vRYXF1erhQEAAMAzlhq6l156SUVFRUpLS3O7hu7MXSMAAABw8Vlq6Hbt2qWPP/5YLVq08FY9AAAA51WX76lqF0vX0F1++eUKCAjwVi0AAACoAUsJ3QMPPKBJkybpvvvu02WXXeY21Nq9e/daLw4AAOBsLFtiZqmhmzt3riQpNTXVbbvL5dL3339fe1UBAADAY5Yaur1793qrDgAAAI9wL1czSw1dZmbmOV9r06bNBRcDAAAA6yw1dFFRUXK5XJWd8S+voWPIFQAAwB6WGrpPP/3U7Xlubq7+/Oc/q3///rVaFAAAwLnU5Ts22MVSQ9e2bVvT80WLFmn48OGKiYmp1cIAAADgGUsN3bnk5+fXxscAAABUi4WFzSw1dCtWrHB7Xlpaqu3btys8PLw2awIAAIAFlm/99Uu+vr7q1q2bJkyYUKtFAQAAwHOWGro333zTW3UAAAB4hDtFmHnU0G3cuLHafYYNG3aBpQAAAKAmPGroli9fft7XXS4XDR0AALgouFOEmUcN3datW71dBwAAAGrI8rIl3333nTZs2KAff/xRLVu21IgRIxQZGemN2gAAAEy4hs7Mx8rOO3bs0OjRo5WXl6errrpKBQUFGj9+vD755BNv1QcAAIBqWEroli9friVLlmjw4MGV27Zs2aKVK1dqwIABtV4cAAAAqmcpofv3v/+tQYMGuW0bNGiQDh06VJs1AQAAnJNh4//qKksNXZMmTbR//363bXv37lXLli1rtSgAAAB4zqMh19LSUvn5+emuu+7SxIkTNWHCBF1xxRU6cuSIXn31VY0ePdrbdQIAAEiSKli2xMSjhu43v/mN7rnnHo0cOVIlJSVKTExUdna22rZtq3vvvVfjx4/3dp0AAAA4B48aumnTpmndunVKTEzU4MGD9dJLL6lr167erg0AAAAe8OgaupEjR+qdd97Rm2++KR8fH40ZM0b33HOPtmzZovLycm/XCAAAUMmw8VFXWZoUER4ermeeeUaff/65BgwYoGXLlql///565ZVXvFUfAAAAqmGpoTujSZMmeuCBB7R+/XrddNNNeuGFF2q7LgAAgCpVyLDtUVfVqKHbtWuXZs6cqb59++qnn37SihUrarsuAAAAeMjjO0Xk5OTonXfe0YYNG3Ts2DHFxMTonXfeUceOHb1ZHwAAgJu6nJTZxaOGbsqUKdq2bZvatGmj0aNHKzY2Vg0bNvR2bQAAAPCARw1dUVGRVqxYob59+8rlcnm7JgAAAFjgUUP3l7/8xdt1AAAAeMTgThEmNZoUAQAAgLrD40kRAAAAdQGTIsxI6AAAAByOhg4AAMDhGHIFAACOYjDkakJCBwAA4HAkdAAAwFFYtsSMhA4AAMDhSOgAAICjsGyJGQkdAACAw9HQAQAAOBxDrgAAwFGYFGFGQgcAAOBwtid0fTOO2l2CLY4V59tdgi2C/QPtLsEWvq76+btTff0tuqjstN0l2CKs02/tLsEWGemb7C6h3mFShFn9/FcGAADgEkJDBwAA4HC2D7kCAABYwb1czUjoAAAAHI6EDgAAOEpFPZ1wdT4kdAAAAA5HQgcAAByFa+jMSOgAAAAcjoYOAADA4RhyBQAAjsKkCDMSOgAAAIcjoQMAAI7CpAgzEjoAAACHo6EDAABwOIZcAQCAozApwoyEDgAAwOFI6AAAgKMwKcKMhA4AAMDhaOgAAAAcjiFXAADgKEyKMCOhAwAAcDgSOgAA4ChMijCzlNBFRUVpxYoV+vHHH71VDwAAACyy1ND94Q9/0HfffadBgwZp3Lhxev/993X69Glv1QYAAGBiGBW2Peoql2FYv7IwJydHSUlJ2rx5sw4fPqzo6Gjdeeeduu666ywX0K6Z9fdcCo4V59tdgi38fevnKL+vq35erlqDHy+XhKKy+vmLbtPAhnaXYIuM9E12l2ALvxYdbDv2lc1vsO3Y/875xrZjn0+N/pVp3ry5hg4dqhEjRqht27Z6++23NWHCBI0YMULff/99bdcIAACA87AUl5SUlOjTTz/Ve++9px07dqhz584aPny4hgwZossuu0wrVqxQXFycPv30U2/VCwAA6rkKJkWYWGrobr75Zvn6+uqOO+7Q+vXr9etf/9rt9ejoaG3cuLE26wMAAEA1LDV0Cxcu1IABA+Tv71/l6506ddLWrVtrpTAAAICq1Nfrc8/HUkMXHR2tjIwMZWVlVX4xS0tLtX//fo0bN84b9QEAAKAalhq6xMREvfDCC3K5XJJ+7pBdLpeuueYaGjoAAACbWGro/va3v2n58uXy9/fX1q1bNWPGDD311FNq3bq1t+oDAABww6QIM0vLluTn52vgwIG6+uqr9d1336lJkyZ68skntXnzZm/VBwAAgGpYSuguv/xyFRQUqFWrVvrhhx9kGIaaNWumEydOeKs+AAAAN0yKMLPU0HXv3l1Tp07VsmXLdO211+r5559XQECAWrVq5a36AAAAUA1LQ66zZs3Sr371K5WVlSk+Pl6ffvqp1q9fr/j4eG/VBwAA4KbCMGx71FU1updrbeJervUL93KtX+rrsAj3cq1fuJfrxde6ybW2Hfs/ef+y7djn49G/ritWrKh2n7i4uAsuBgAAANZ51NDt2rVLknTq1Cl9++23uvbaa3XFFVcoKytL33zzjW655RavFgkAAHCGwbIlJh41dG+++aakn6+hGzFihO65557K195991394x//8E51AAAAqJalC3v+8Y9/aNSoUW7bYmJitHPnzlotCgAA4FwMw7DtUVdZauiaNWum1NRUt207duzQ5ZdfXqtFAQAAwHOWphxOmDBBDz30kAYNGqQ2bdooIyNDn3zyiZYsWeKt+gAAAFANSw3dXXfdpSuuuEJJSUn67rvvFBoaqlWrVunGG2/0Vn0AAABuuJermeVFwXr16qUuXbroxx9/VMuWLdW6dWtv1AUAAOBoOTk5mjNnjlJSUuTr66uYmBg9/vjjatDA3H59/vnnWrp0qTIyMtS6dWv94Q9/UL9+/Tw+lqVr6AoKCjR58mT17t1bo0aNUlRUlH73u98pP79+LpILAAAuPqdMipg+fbqCg4O1fft2bdiwQV988YVWrVpl2u/QoUOaMmWKpk2bpi+//FJTpkzR9OnTlZWV5fGxLDV0zz33nAoLC/XBBx/o66+/1nvvvaeKigo9++yzVj4GAADgknb48GGlpKToscceU1BQkMLCwjRp0iStWbPGtO+7776ryMhIDRgwQA0aNFB0dLS6d++u9evXe3w8Sw3dtm3b9Nxzz6ljx44KCAhQly5d9Oyzz+qTTz6x8jEAAAA1Zue9XEtKSlRQUOD2KCkpMdV44MABNWnSRK1atarc1rFjR2VmZppGNtPT09WlSxe3bZ06ddLevXs9/ppYauiKi4vVqFEjt22NGzdWRUWFlY8BAABwpMTEREVERLg9EhMTTfsVFhYqKCjIbduZ50VFRdXuGxgYaNrvfDyaFLF7925FRETohhtu0Isvvqjf//73crlcMgxDL774oq677jqPDwgAAOBUEyZM0Pjx4922+fv7m/YLDg5WcXGx27Yzz0NCQty2BwUF6dSpU27bTp06ZdrvfDxq6B566CHt2bNHM2fO1H333aekpCS1bdtWP/74o1wul/761796fEAAAIALYecdG/z9/ats4M7WuXNn5eXlKTs7Wy1atJAkHTx4UKGhoabRzi5duigtLc1tW3p6urp27epxXR4NuZ75wl111VX66KOPNG3aNEVFRWnWrFnavHmzOnbs6PEBAQAALnXt27dXRESEFi9erIKCAmVkZGjlypWKjY017RsTE6OUlBRt3rxZZWVl2rx5s1JSUjR06FCPj+cyPGhzb7zxRu3Zs8famXioXbP6OVx7rLh+LvXi72t56cNLgq/L0uWql4y6fN9DbyoqO213CbZoGtjQ7hJskZG+ye4SbOHXooNtx76soX1B0omCgx7vm52drYULF2rXrl3y8fHRsGHD9Pvf/16+vr7q1q2bFixYoJiYGEnS9u3btXTpUh05ckRt27bVY489pr59+3p8LI8aumuuuUZt2rQ57z6ffvqpxwf9JRq6+oWGrn6hoatfaOjqFxq6usWjf139/PwUFxfn7VoAAABQAx41dA0aNNDw4cO9XQsAAEC16mv6fz6WJkUAAACg7vEooTtzwR4AAIDdKgiaTDxK6BYsWODtOgAAAFBD9XPKIQAAcCxDJHRnq59rKQAAAFxCaOgAAAAcjiFXAADgKEyKMCOhAwAAcDgSOgAA4Cisj2tGQgcAAOBwNHQAAAAOx5ArAABwFNahMyOhAwAAcDgSOgAA4ChMijAjoQMAAHA4EjoAAOAoJHRmJHQAAAAOR0MHAADgcAy5AgAAR2HA1YyEDgAAwOFcBlcWAgAAOBoJHQAAgMPR0AEAADgcDR0AAIDD0dABAAA4HA0dAACAw9HQAQAAOBwNHQAAgMPR0AEAADgcDR0AAIDDXRIN3fz583XLLbcoJyfHbXtZWZlGjhypCRMm6FK7IcZVV12lXbt2efUYxcXFGjVqlN555x2vHscKb573Dz/8oLi4ON10003q2bOnJk2apIyMDK8cyypvnvf333+v++67TxEREerZs6cee+wxHT9+3CvHsupi/D2XpMcee0xjx471+nE85c3z/uabb3T11VerW7dulY8xY8Z45Vhn+93vfqe4uLgqX3vrrbd08803q6Sk5Jzvj4qKqlM/j2pq7ty5lV/76667zvT9+PLLLz3+rLFjx+qll17yYrVwikuioXviiSfUokULPfHEE27bX3rpJWVnZ2vJkiVyuVw2VedMBw4c0JgxY/T111/bXcpFM3nyZF122WXaunWrtm7dqiZNmmjSpEl2l+VVJSUleuihh9SzZ0/t2rVLH3/8sY4dO6ZnnnnG7tIumg0bNuiDDz6wu4yL5ttvv1X37t311VdfVT7WrFlzUY49duxYbdu2TceOHTO9tnbtWt19993y9/e/KLXYaeHChZVf+wULFqhNmzZu34/IyEi7S4QDXRINXUBAgF544QWlpqbqzTfflCSlpKRo1apVWrZsmfLz8/XII4+oZ8+e6tevn1544YXK3wINw9Arr7yiIUOGKDIyUt27d9fMmTN16tQpSdKsWbM0depUDR48WDfddJOOHDli23l6qqSkREuWLNHgwYPVrVs39erVS0899VRlSnnq1CnNmzdPPXr0UN++fbVs2TJFRUVVJgJffPGF7r//fg0fPlxt2rSx81QsuZDzPnHihFq0aKFp06YpODhYISEhuu+++7R//36dOHHC5jM7vws5b39/f/3jH//QxIkT1aBBA504cULFxcVq1qyZzWdVvQv9ey5J6enpWrlype666y67TsOyCz3vb7/9Vl27drWl9r59+6pNmzZ699133bZ//fXXOnDggEaNGqVXXnlFAwYMUEREhGJjY7V9+/YqP+vsZOqHH37QVVddpR9++EHSzynn+vXrNWjQIN1www165JFH9N133+nuu+9Wt27ddOedd+rw4cOV79+0aZOGDBmiiIgIjRgxQjt27PDCV6B6+/bt00MPPaQePXro1ltv1fz583Xy5MnK1//+97+rf//+6tatmx5//HEVFxdXvlZQUKDZs2dr4MCBCg8PV58+ffSnP/1J0s/nFxERodOnT1fu/+GHH6pfv36X3AhWvWVcQt555x3jhhtuMP71r38Z/fr1M1avXm0UFhYa/fr1M5YuXWqcOnXKyMzMNGJjY42lS5cahmEYmzZtMm655Rbj3//+t2EYhpGenm706NHDeOuttwzDMIzHH3/cCA8PN/bt22ecOHHCrlMz6dKli7Fz584qX3vllVeM3/72t0ZWVpZhGIaxZ88e49prrzWSk5MNwzCMOXPmGMOHDzcyMzONgoIC47HHHnP7vNzcXOPUqVOGYRhGv379jLfffvsinJFnvHneZ1u2bJnRr18/75yIRRfjvEeNGmV06dLFiI6ONo4dO+bdE/KQN8+7uLjYuOOOO4zPPvvMWL58uXHvvfdenJPygDfPe/DgwcbYsWON2267zejVq5cxbdo04z//+c/FOTHDMF577TVjwIABRkVFReW2P/zhD8bMmTON5cuXG7feeqvx3XffGaWlpcamTZuMrl27Gt98841hGO4/j+69915j+fLllZ+RkZFhdOnSxcjIyDAM4+ev4ZgxY4zjx48bWVlZRmRkpNGnTx8jPT3dKCwsNO6++25j1qxZhmEYxmeffWZEREQYKSkpRllZmbF161YjPDzc2L9/v9e/Hm+//Xblz5nc3FyjR48exjPPPGMUFxcbP/30k3HfffcZjzzyiGEYhpGcnGx07drVSE5ONkpLS401a9YYXbp0qfw6zJs3z7j//vuNEydOGBUVFcaHH35odOnSxTh06JBx+vRpo3v37samTZsqjz1hwgTjxRdf9Po54uK4JBK6M4YPH66BAwdW/gY2ZswYffbZZyopKdGMGTMUEBCg1q1ba9q0aZVDDLfeeqs2bNig9u3bKzc3V8ePH1eTJk2UlZVV+bnh4eHq0qWLGjdubNepWTJy5EitWrVKLVu21E8//aRTp04pJCREWVlZKi0tVVJSkh599FG1bt1aISEhmjt3rnx9fSvf37RpUwUEBNh4BjVzoef9S2vXrtVrr72mRYsWXeSzsK62znvVqlVKSUlRly5dNH78eJWXl9twNp670PNeuHChbrnlFvXt29fGs7DuQs67vLxcl19+uXr37q23335bH3zwgVwulx5++OGL9v2OjY1Vdna2du7cKUnKy8vTli1bdN999+ntt9/Www8/rF//+tdq0KCBoqOjFRUVpQ0bNtToWPfee6+aNGmiyy+/XJ07d9bAgQPVsWNHBQcH66abbtKPP/4oSVq9erXuuecede/eXb6+vurXr5+ioqK0bt26WjtvT3z66afy8/PT73//ewUGBqply5aaM2eOtm7dqmPHjikpKUkDBw5Ur1691KBBA40ePVrXXntt5funTJmiZcuWqWHDhjp69Gjlz/GffvpJ/v7+uuOOO/Tee+9JknJycrRjxw4NHz78op4jvKeB3QXUtri4OL333nuaNm2aJOnHH39Ubm6uunfvXrmPYRgqLS1VTk6O/P399cILL2jbtm1q1qyZrrnmGpWWlrpF0JdffvlFP48LUVxcrIULFyo1NVWhoaG69tprZRiGKioqlJeXp+LiYrVt27Zy/4YNG6pp06Y2Vlw7auO8S0pK9PTTT2vz5s1KTEzUTTfddLFPw7La+n4HBgYqMDBQs2fP1s0336x9+/a5/WNR11zIeSclJWnv3r0X/R/s2nAh5+3r66tVq1a5fd6cOXPUq1cvHTx4UF26dPF6/Y0aNVJMTIz+/ve/q1evXnr77bd17bXX6vrrr1d2drbCwsLc9r/iiiu0d+/eGh2rSZMmlX/29fXVZZddVvncx8en8uf8jz/+qJSUFK1du7by9fLy8ov+339OTo7atGnj9ovHFVdcUVljVlaWfv3rX7u955dfr5ycHCUkJOhf//qXrrjiisqh9YqKCknSiBEjNGrUKOXk5CgpKUk33nij6esN57rkGjofHx+3/w8NDVW7du304YcfVu5TUFCgnJwcNWvWTPPnz1dmZqa2bt2qhg0bSpKGDBni9plOm1Axe/ZsXXbZZdqxY4cCAgJUUVFR2dA2b95cgYGByszMVIcOHSRJRUVFdWZW44W40PPOzc3VxIkTVVJSog0bNjjmB92FnPcPP/yg++67T+vWrav8xeXM9aW//MevLrqQ837vvff073//WzfffLMk6fTp0yovL1dkZKSSkpLq9LWjF3Le//nPf7Rq1SpNnTpVISEhkv7v+x0YGHjRzmHs2LEaPny4jh8/rrfeektTp06VJLVt29Y0szwjI6PKX6p9fHxUWlpa+byqn2Ge/uwODQ3VsGHD9PDDD1duy8zMvKhfE+nn88/MzFR5eXllU3fmuu2WLVsqNDTU9PU5evSoOnfuLEmaNm2aoqKi9Je//EUNGjSo/Pqe0bVrV3Xq1EkfffSRNm3aVKdmduPCXVJDrlXp16+fCgsL9ec//1klJSXKz8/X448/rkcffVQul0sFBQUKCAiQr6+vTp8+rddee0379+93+0FRV+Xm5uro0aNuj7Kysspz8vHxUUFBgf7rv/5LBQUFKi0tlY+Pj2JjY/XSSy8pKytLxcXFevrpp+v88NoveeO8S0tL9eCDD6phw4Zau3ZtnWzmvHHebdu2VZMmTfT000+rsLBQubm5WrBggW699Va3lMdO3jjvv/zlL/rqq6/05Zdf6ssvv9TDDz+siIgIffnll3WmmfPGeTdt2lSbNm3SCy+8oNOnT1d+v3v16qV27dpdtHPr1KmTIiIi9Mwzz6i4uFgDBw6UJN1111165ZVXlJaWpvLycm3ZskVbt26tcliwY8eO2r59u/Lz83Xy5Em9+uqrNa5n5MiReuONN/TPf/5T0s8TR0aMGHHRZz+fGf5funSpTp06pWPHjikhIUE33XST2rZtqzvvvFOffPKJtm3bprKyMr377rv65ptvKt9/8uRJBQYGytfXV7m5uZWXjPzy37MRI0borbfe0qFDhyq/7rg0XPINXcOGDbVq1Srt2rVLt956qwYMGCAfHx+9/PLLkqTp06fr1KlTuvnmmxUVFaWvv/5aQ4cO1f79+22uvHrTp09X37593R6HDx/W7NmztXfvXvXo0UO33367CgoK1KdPn8pzmjlzpjp06KDo6GgNGjRIoaGh8vHxkZ+fn81n5BlvnPe2bduUlpam1NRU9erVy21NqMzMTJvP+GfeOG+Xy6WVK1eqrKxMUVFRGjp0qFq3bq3nn3/e5rP9P/w9r73zDgwM1J///GcdPHhQvXv31qBBg9SwYUMtW7bsop/fvffeq40bN+qee+6p/J6MHz9eY8aM0aOPPqrIyEglJibq+eefV48ePUzvnzBhgpo3b67+/ftr6NChioqKqnEtt99+u2bMmKH4+HjdeOONmjZtmsaNG3fRE6xGjRrpr3/9q/bv36++ffvqjjvuUNu2bfXiiy9KkiIiIvRf//VfeuaZZxQZGamPPvpIt9xyS+X7z1wucuONN2rEiBFq1aqVrr32Wrd/z4YMGaL09HRFR0crKCjoop4fvMtlGMxXrm9SU1N11VVXVU7yKCgoUEREhD766CO1b9/e3uK8iPPmvDnv9vYWB9uVl5erd+/e+tOf/qQbbrjB7nJQiy65a+hQvddee02NGzfWggUL5HK5tHz5cl155ZWX/A97zpvz5rxRnx04cEBbtmxRaGgozdwliISuHsrKytKCBQu0e/dulZeXKyIiQk8++eRFvYbGDpw35815oz47Myy9fPly2xaXhvfQ0AEAADjcJT8pAgAA4FJHQwcAAOBwNHQAAAAOR0MHAADgcDR0AAAADkdDBwAA4HA0dAAAAA5HQwcAAOBwNHQAAAAO9/8Ap6iZN/XGvakAAAAASUVORK5CYII=\n"
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# the best way to visualize corerelations matrices is heatmap\n",
    "plt.figure(figsize = (8,8))\n",
    "sns.heatmap(corr)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### We can see from the above figure that correlations between veriablesa excpet Volume and Year are close to zero."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-07-05T12:00:37.606037Z",
     "start_time": "2023-07-05T12:00:37.501298Z"
    }
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "          Volume     Today\n",
      "Year    0.539006  0.030095\n",
      "Lag1    0.040910 -0.026155\n",
      "Lag2   -0.043383 -0.010250\n",
      "Lag3   -0.041824 -0.002448\n",
      "Lag4   -0.048414 -0.006900\n",
      "Lag5   -0.022002 -0.034860\n",
      "Volume  1.000000  0.014592\n",
      "Today   0.014592  1.000000\n"
     ]
    }
   ],
   "source": [
    "print(corr.iloc[:,-2:])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "ExecuteTime": {
     "end_time": "2023-07-05T12:00:37.889411Z",
     "start_time": "2023-07-05T12:00:37.512370Z"
    }
   },
   "outputs": [
    {
     "ename": "TypeError",
     "evalue": "scatterplot() takes from 0 to 1 positional arguments but 2 were given",
     "output_type": "error",
     "traceback": [
      "\u001B[0;31m---------------------------------------------------------------------------\u001B[0m",
      "\u001B[0;31mTypeError\u001B[0m                                 Traceback (most recent call last)",
      "Cell \u001B[0;32mIn[8], line 3\u001B[0m\n\u001B[1;32m      1\u001B[0m \u001B[38;5;66;03m# we choose the x axis as index, chossing year will give a discrete plot\u001B[39;00m\n\u001B[1;32m      2\u001B[0m plt\u001B[38;5;241m.\u001B[39mfigure(figsize \u001B[38;5;241m=\u001B[39m (\u001B[38;5;241m12\u001B[39m,\u001B[38;5;241m6\u001B[39m))\n\u001B[0;32m----> 3\u001B[0m \u001B[43msns\u001B[49m\u001B[38;5;241;43m.\u001B[39;49m\u001B[43mscatterplot\u001B[49m\u001B[43m(\u001B[49m\u001B[43mdata\u001B[49m\u001B[38;5;241;43m.\u001B[39;49m\u001B[43mindex\u001B[49m\u001B[43m,\u001B[49m\u001B[43mdata\u001B[49m\u001B[43m[\u001B[49m\u001B[38;5;124;43m'\u001B[39;49m\u001B[38;5;124;43mVolume\u001B[39;49m\u001B[38;5;124;43m'\u001B[39;49m\u001B[43m]\u001B[49m\u001B[43m)\u001B[49m\n",
      "\u001B[0;31mTypeError\u001B[0m: scatterplot() takes from 0 to 1 positional arguments but 2 were given"
     ]
    },
    {
     "data": {
      "text/plain": "<Figure size 1200x600 with 0 Axes>"
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# we choose the x axis as index, chossing year will give a discrete plot\n",
    "plt.figure(figsize = (12,6))\n",
    "sns.scatterplot(data.index,data['Volume'])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4.6.2 Logistic Regression"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#store a copy of the original data, as we are going to make some changes in the data\n",
    "data_orig = data.copy()\n",
    "data['Direction'] = data['Direction'].map({'Down':0,'Up':1})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#logit = sm.Logit('Direction ~ Lag1 + Lag2 + Lag3 + Lag4 + Lag5 + Volume',data = data)\n",
    "#the above line was giving error\n",
    "\n",
    "X = data.iloc[:,1:-2]\n",
    "X = sm.add_constant(X)\n",
    "y = data['Direction']\n",
    "results = sm.Logit(y,X).fit()\n",
    "print(results.summary())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# probablitie that the given oberstaion will result in 'UP'\n",
    "prob = results.predict(X)[:10]\n",
    "print(list(prob))\n",
    "\n",
    "predicted_classes = np.where(prob <=0.5,'Down','Up')\n",
    "pd.DataFrame({'Probabilities':prob,'Classes':predicted_classes})"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# confusion matrix and accuracy\n",
    "# Using sklearn\n",
    "lr = LogisticRegression()\n",
    "lr.fit(X,y)\n",
    "lr.coef_\n",
    "# we can see that sklearn and statsmodel, don't have the exact same coefficients"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cm = confusion_matrix(y,lr.predict(X))\n",
    "accuracy = accuracy_score(y,lr.predict(X))\n",
    "\n",
    "print('Accuracy is ',accuracy)\n",
    "print('confusion_matrix - ')\n",
    "print(cm)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "NOTE - Confusion matrix drawn in the book has actual values as columns, and predicted values as rows.\n",
    "       While in Sklearn, actual values are as rows, and predicted values are as columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#using statsmodel\n",
    "table = results.pred_table(threshold=0.5)\n",
    "table"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Computin train and test accuracy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "train = data[data['Year']<2005]\n",
    "test = data[data['Year'] == 2005]\n",
    "print('Shape of train is',train.shape)\n",
    "print('Shape of test is ',test.shape)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_train = train.iloc[:,1:-2]\n",
    "X_train = sm.add_constant(X_train)\n",
    "y_train = train['Direction']\n",
    "\n",
    "X_test = test.iloc[:,1:-2]\n",
    "X_test = sm.add_constant(X_test)\n",
    "y_test = test['Direction']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "results = sm.Logit(y_train,X_train).fit()\n",
    "pred_test = results.predict(X_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pred_classes = np.where(pred_test > 0.5,1,0)\n",
    "cm = confusion_matrix(y_test,pred_classes)\n",
    "cm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "accuracy = accuracy_score(y_test,pred_classes)\n",
    "print('Accuracy is ',accuracy)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# The accuracy is just 48%, even random guessing will do a better job."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Using only Lag1 and lag 2\n",
    "predictors = ['Lag1','Lag2']\n",
    "X_train = train[predictors]\n",
    "X_train = sm.add_constant(X_train)\n",
    "\n",
    "X_test = test[predictors]\n",
    "X_test = sm.add_constant(X_test)\n",
    "# y_train and y_test remains the same"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "results = sm.Logit(y_train,X_train).fit()\n",
    "pred_test = results.predict(X_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pred_classes = np.where(pred_test > 0.5,1,0)\n",
    "cm = confusion_matrix(y_test,pred_classes)\n",
    "cm"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "accuracy = accuracy_score(y_test,pred_classes)\n",
    "print('Accuracy is ',accuracy)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# the accuracy is better that before"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# prediction\n",
    "results.predict([[1,1.2,1.1],\n",
    "                 [1,1.5,-0.8]])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4.6.3 Linear Discriment Analysis"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# We use sklearn for LDA\n",
    "#from sklearn.discriminant_analysis import LinearDiscriminantAnalysis as LDA"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "https://scikit-learn.org/stable/modules/generated/sklearn.discriminant_analysis.LinearDiscriminantAnalysis.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "predictors = ['Lag1','Lag2']\n",
    "\n",
    "X_train = train[predictors]\n",
    "X_test = test[predictors]\n",
    "\n",
    "lda = LDA()\n",
    "# in the parameters we can also provide prior probabilities represented by pie in chapter\n",
    "lda.fit(X_train,y_train)\n",
    "pred = lda.predict(X_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Prior Probs are - ',lda.priors_)\n",
    "print('Class Means are  - ',lda.means_)\n",
    "print('Coeff are - ',lda.coef_)\n",
    "\n",
    "# Mean is the average value of observations for each class"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "cm = confusion_matrix(y_test,pred)\n",
    "print(cm)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Accuracy using LDA is ',accuracy_score(y_test,pred))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "plt.figure(figsize = (12,8))\n",
    "tmp = pd.DataFrame({'Lag1':X_train['Lag1'],'Lag2':X_train['Lag2'],'Direction':y_train})\n",
    "sns.scatterplot(y = 'Lag2',x = 'Lag1',hue = 'Direction',data = tmp)\n",
    "\n",
    "x1 = X_train['Lag1']\n",
    "b, w1, w2 = lda.intercept_[0], lda.coef_[0][0], lda.coef_[0][1]\n",
    "y1 = -(b+x1*w1)/w2    \n",
    "plt.plot(x1,y1)\n",
    "\n",
    "#mean \n",
    "plt.scatter(lda.means_[0][0],lda.means_[0][1],marker = \"+\",color = 'blue',s = 800)\n",
    "plt.scatter(lda.means_[1][0],lda.means_[1][1],marker = \"+\",color = 'red',s = 800)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4.6.4 Quadratic Discriminant Analysis"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "https://scikit-learn.org/stable/modules/generated/sklearn.discriminant_analysis.QuadraticDiscriminantAnalysis.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "qda = QDA()\n",
    "qda.fit(X_train,y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Mean for class 0 is - ',qda.means_[0])\n",
    "print('Mean for class 1 is - ',qda.means_[1])\n",
    "print('Prior probalbilities - ',qda.priors_)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pred = qda.predict(X_test)\n",
    "cm = confusion_matrix(y_test,pred)\n",
    "print(cm)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Accuracy using QDA is ',accuracy_score(y_test,pred))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# QDA has the best performance with almost 60% accuracy\n",
    "# We can see from the graph that there is no possible linear spearation, and non linear models will perform better"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4.6.5 K-Nearest Neighbors"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "knn = KNeighborsClassifier(n_neighbors=1)\n",
    "knn.fit(X_train,y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pred = knn.predict(X_test)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print('Accuracy using KNN for 1 nieghbor is ',accuracy_score(y_test,pred))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "knn = KNeighborsClassifier(n_neighbors=3)\n",
    "knn.fit(X_train,y_train)\n",
    "pred = knn.predict(X_test)\n",
    "print('Accuracy using KNN for 1 nieghbor is ',accuracy_score(y_test,pred))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Even if we increase n to 3, accuracy is 53.17"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 4.6.6 An Application to Caravan Insurance Data\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data = pd.read_csv('E:\\programming\\dataset\\Into_to_statstical_learning\\Caravan.csv',index_col=0)\n",
    "print(data.shape)\n",
    "data.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data['Purchase'].value_counts() / len(data)\n",
    "# Only around 6% of the people purchased the insurance."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# From a sample of means if the observations it can be seen that the variables falls in different ranges\n",
    "# We want all the variables to have a same range. \n",
    "data.mean()[:5]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# we will applt standard scaler which will transform the data such that it has 0  mean, and std deviation of 1\n",
    "# but before that lets split the data into test and train data "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "test = data.iloc[:1000,:]\n",
    "train = data.iloc[1000:,:]\n",
    "\n",
    "X_train = train.drop('Purchase',axis = 1)\n",
    "y_train = train['Purchase']\n",
    "X_test = test.drop('Purchase',axis = 1)\n",
    "y_test = test['Purchase']"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Applying the standard Scaler\n",
    "scaler = StandardScaler()\n",
    "\n",
    "X_train = scaler.fit_transform(X_train)\n",
    "X_test = scaler.transform(X_test)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "https://scikit-learn.org/stable/modules/generated/sklearn.preprocessing.StandardScaler.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# changing labels from string to binary\n",
    "y_train = np.where(y_train == 'Yes',1,0)\n",
    "y_test = np.where(y_test == 'Yes',1,0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# KNN with n_neighbors = 1\n",
    "knn_1 = KNeighborsClassifier(n_neighbors=1)\n",
    "\n",
    "knn_1.fit(X_train,y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "pred = knn_1.predict(X_test)\n",
    "cm = confusion_matrix(y_test,pred)\n",
    "print('The accuracy for k = 1 is ',accuracy_score(y_test,pred))\n",
    "print('Confusion Matrix')\n",
    "pd.DataFrame(cm,columns = ['No','Yes'],index = ['No','Yes'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Out of 76 predicted values, 9 are in True in real.\n",
    "# accuracy = 11%, which is better than random guesing"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# KNN with n_neighbors = 3\n",
    "knn_3 = KNeighborsClassifier(n_neighbors=3)\n",
    "\n",
    "knn_3.fit(X_train,y_train)\n",
    "\n",
    "pred = knn_3.predict(X_test)\n",
    "cm = confusion_matrix(y_test,pred)\n",
    "print('The accuracy for k = 3 is ',accuracy_score(y_test,pred))\n",
    "print('Confusion Matrix')\n",
    "pd.DataFrame(cm,columns = ['No','Yes'],index = ['No','Yes'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# for 27 totol predicted values, 6 are Yes in real\n",
    "# Accuracy is 22%"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# KNN with n_neighbors = 5\n",
    "knn_5 = KNeighborsClassifier(n_neighbors=5)\n",
    "\n",
    "knn_5.fit(X_train,y_train)\n",
    "\n",
    "pred = knn_5.predict(X_test)\n",
    "cm = confusion_matrix(y_test,pred)\n",
    "print('The accuracy for k = 5 is ',accuracy_score(y_test,pred))\n",
    "print('Confusion Matrix')\n",
    "pd.DataFrame(cm,columns = ['No','Yes'],index = ['No','Yes'])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# for k = 5, accuracy is 33.3%"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Additional  - Finding the optimal value of k"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def run_knn(n_neighbors):\n",
    "    accuracy_list = []\n",
    "    for n in n_neighbors:\n",
    "        knn = KNeighborsClassifier(n_neighbors=n)\n",
    "\n",
    "        knn.fit(X_train,y_train)\n",
    "\n",
    "        pred = knn.predict(X_test)\n",
    "        accuracy_list.append(accuracy_score(y_test,pred))\n",
    "    return accuracy_list    "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "neighbors = [i for i in range(1,10)]\n",
    "accuracy_list = run_knn(neighbors)\n",
    "plt.plot(neighbors,accuracy_list)\n",
    "plt.xlabel('N_neighbors')\n",
    "plt.ylabel('Test Accuracy')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
